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ABSTRACT 

The onset of convection induced by coupled surface 
tension gradient and buoyancy forces is investigated with 
temperature dependent viscosity. Both surface tension and 
viscosity are assumed to vary linearly with temperature. The 
linear viscosity approximation was found to be an adequate 
representation of relevant published experimental data. One 
limiting case, Ma = 0, is the buoyancy driven convection 
problem typically referred to as the Rayleigh-Benard stability 
problem. The other limiting case, Ra = 0, is the surface tension 
(gradient) driven flow problem referred to as the Marangoni- 
Benard problem. Several investigations of the Rayleigh-Benard 
problem with temperature dependent viscosity have been 
reported. One stability investigation of the Marangoni-Benard 
problem with a linear viscosity variation has also been 
published. Results from that Marangoni-Benard study indicate 
that linear viscosity variation is destabilizing while similar 
studies for the variable-viscosity-Rayleigh-Benard problem 
suggest that the above viscosity profile is stabilizing. In this 
study the variable viscosity analysis is extended to the coupled 
problem, which bridges the above limiting cases. 

The equations and boundary conditions obtained from the 
linear analysis are solved numerically as a generalized 
eigenvalue problem. Neutral stability curves for different 
viscosity slopes are generated for the Marangoni-Benard and 
Rayleigh-Benard problems. It is shown that the curves can be 
collapsed to a single curve by appropriately scaling the results 
for each of the limiting cases. The critical Marangoni number is 
determined as a function of the slope of the viscosity 
temperature variation, e, for different values of the Rayleigh 
number. Regression analyses of the numerical results are 
performed to provide a convenient means of computing the 
critical Marangoni number as a function of e and the Rayleigh 
number. The linear temperature-dependent viscosity considered 
in this analysis gives results which are consistent with the 
Rayleigh-Benard studies employing an exponentially 
temperature-dependent viscosity. When the viscosity decreases 
linearly with temperature, the coupled buoyancy-surface 
tension problem, including the limiting special cases of Ra = 0 
and Ma = 0, is found to be more stable than the constant 


viscosity case. The difference between this and the previous 
Marangoni-Benard study is explained. 

NOMENCLATURE 

a-f finite difference coefficients (Eqs. (21) and (22)) 

Bi s surface Biot number 

c v c 2 proportionality constants for AT C 

d fluid layer thickness or depth 

g gravitational acceleration (or finite diff. coef for 

Eqs. (21) and (22) 

h surface heat flux temperature gradient 

i 

k thermal conductivity 

k unit vector in z direction 

Ma Marangoni number 

P pressure 

p pressure disturbance 

Pr Prandtl number 

Q surface heat flux 

Ra Rayleigh number 

T temperature 

t time 

AT temperature difference 

U velocity vector 

w velocity disturbance amplitude 

u velocity disturbance vector 

x,y,z cartesian coordinates 

Greek symbols 
a wavenumber 

(3 base temperature gradient 

e viscosity - temperature slope 

k thermal diffusivity 

X eigenvalue (dimensionless frequency) 

p dynamic viscosity 

v kinematic viscosity 

0 temperature disturbance 

<j) temperature disturbance amplitude 



p fluid density 

a surface tension 

£ volumetric thermal expansion coefficient 

Subscripts 

0 parameter/variable value at lower surface 

b base quantity 

CL result based on Cloot and Lebon '85 

c critical value 

d parameter/variable value across fluid layer 

max maximum value 

min minimu m value 

s property/parameter value at upper surface 

Superscripts 

mb Marangoni-Benard (Ra=0) 

rb Rayleigh-Benard (Ma=0) 

T Transpose (see eqs. 2 and 11) 

* dimensional quantity 

INTRODUCTION 

The onset of cellular convection with temperature 
gradients imposed normal to the free surface continues to be an 
active research topic as illustrated by Koschmeider’s (1993) 
extensive literature survey. Microgravity science applications 
of thermocapillary flows have motivated the study of surface 
tension induced instability known as the Marangoni-Benard 
problem in both one-g and microgravity environments 
(Legros, J.C., et al. 1990). In space, the effect of buoyancy is 
minimal and one can consider surface tension to be the sole 
driving force. Then the temperature difference necessary to 
achieve the critical Marangoni number, Ma c , depends only on 
the fluid properties and fluid thickness. On the other hand, in 
ground-based studies of surface tension driven problem, Ma c 
is also influenced by buoyancy, that is by the Rayleigh number. 
To minimize the effect of buoyancy in the one-g environment, 
the liquid layer must be thin but this means that a larger 
temperature difference must be imposed across the layer of 
fluid in order to reach the critical Marangoni number. As the 
temperature difference increases it is no longer appropriate to 
use constant values for the fluid properties. The effect of 
variable fluid properties on Ma c becomes an important 
consideration. 

The strong dependence of viscosity on temperature for 
most fluids, has prompted several studies of the onset of 
cellular convection which include this feature, (Palm, 1960, 
Stengel et. al.,1982. White, 1988, and Koschmeider, 1993). 
The studies almost exclusively consider the buoyancy induced 
instability problem. In fact, only one study treats the surface 
tension driven instability with a temperature dependent 
viscosity (Cloot and Lebon, 1985). For small to moderate 
viscosity variations, Stengel et al. and White found that the 
variable viscosity had a stabilizing effect on the Rayleigh- 
Benard problem while Cloot and Lebon concluded that the 


viscosity variation is destabilizing for the Marangoni-Benard 
problem. In this paper we investigate the effects of variable 
viscosity on the combined buoyancy and surface tension driven 
instability which was originally considered by Nield (1964). 
The Marangoni-Benard and Rayleigh-Benard instabilities are 
discussed as special cases. Viscosity is assumed to depend 
linearly on temperature throughout the study. 

The governing equations and boundary conditions are 
developed in the following section. The numerical solution of 
the resulting eigenvalue problem is briefly discussed. The 
neutral stability curves for several viscosity-temperature slopes 
are presented and these show the stabilizing effect of variable 
viscosity for both the Marangoni-Benard and Rayleigh-Benard 
instabilities (for small viscosity variations). After rescaling the 
neutral stability results with the limiting critical parameter 
values, the variable viscosity curves are shown to collapse onto 
the constant viscosity curve. Variable viscosity results for the 
coupled surface tension and buoyancy induced instability are 
then compared with Nield’s constant viscosity results. The 
relative contribution of the surface tension and buoyancy in 
establishing the instability is shown as a function of fluid depth 
for a high Prandtl number fluid, silicone oil, and a low Prandtl 
number fluid, mercury. The critical temperature difference 
across the fluid layer indicating that the fluid has lost linear 
stability, AT C , and the corresponding viscosity-temperature 
slope, e, are shown as functions of the fluid depth for both 
fluids. The validity of the linear viscosity profile which 
depends on AT C is then examined as a function of the layer 
thickness, d. 

DEVELOPMENT OF EQUATIONS 

A temperature difference is imposed normal to the free 
surface of a thin liquid layer of fluid of infinite horizontal 
extent and finite thickness, d, as shown in Figure 1. The initial 
steady state or base state of the system is one of no fluid 
motion, with a linear temperature profile across the layer. The 
velocity and temperature profiles illustrated in Fig. 1 can 
immediately be expressed as, U£ = 0 and Tb - T b0 - pz . 
Using the notation of Pearson (1958) and Chandrasehkar 
(1981), U b and T b are respectively, the base flow velocity and 
temperature where the asterisk ’**" denotes dimensional 
quantities and the subscript ”b” denotes that this solution is the 
base flow. The temperature gradient of the base state, P is 
defined as P = -dT b /dz*or P = AT d * /d where 

AT d *=T b ’ 0 -T b V The lower surface is rigid and is held at a 
constant temperature. The upper surface is free and exchanges 
heat with the environment . It is well known that the 
deformation of the free surface also affects the critical 
temperature difference which leads to fluid motion, but in this 
study we focus on the variable viscosity effect and assume that 
the free surface is flat First we give the nondimensional form 
of the governing equations and in the next section we linearize 
about the base state just described in order to determine 
whether small disturbances to the base state will grow or decay. 


2 



Specifically we are interested in the critical values of the 
nondimensional parameters where the change of stability 
occurs. 



Figure 1 . Base State For Thin Liquid Layer Of Infinite Extent 


Nondimensional forms of mass, momentum, and energy 
equations for an incompressible fluid with the Boussinesq 
approximation are given in Eqs. (1) to (3). The derivation of 
these equations with the Boussinesq approximation and 
constant viscosity and their subsequent nondimensionalization 
are well known and we refer the interested reader to 
Chandrasekar (1981), and Drazin & Reid (1982) for details. 
Since dynamic viscosity is the fluid property which depends 
most strongly on temperature, in our study we permit the 
dynamic viscosity to vary with temperature while all other 
thermophysical properties apart from density and surface 
tension, are assumed constant. The inclusion of variable 
viscosity into the governing equations is discussed by Stengel 
et. al. (1982). 


VU = 0 ( 1 ) 

^jU-VP+k RaPr(T-T b0 )+ 

( 2 ) 

PrV (n(T>rVU+(vu) T J 

— = V 2 T (3) 

Dt 

0 , T, P, t are the velocity vector, temperature, pressure, 
and time respectively. The reference values used to 
nondimensionalize the variables; length, velocity, temperature, 
pressure, and time are d, Kj/d, (3d, pjKo/d 2 , d 2 /^, 
respectively. Po is the fluid density and Kq is the fluid thermal 
diffusivity. The subscript 0 indicates that the properties are 
chosen at the lower surface temperature, ^bO. The characteristic 
value of the dynamic viscosity of the fluid, (X, is denoted as 
These reference values are consistent with those used in 
the buoyancy instability studies presented in Chandrasekhar 
(1981) and Drazin and Reid (1982), and the surface tension 
instability investigations of Pearson (1958) and Scriven and 
Stemling (1964). Two dimensionless groups appear in the 
momentum equation, the Prandtl number, Pr, and the Rayleigh 
number, Ra, which are defined as follows: 


Pr = 


Ho 

Po K o 


Ra = 


PoPd 4 ^g 

*0^0 


is the volumetric thermal expansion coefficient and g is 
gravitational acceleration in the negative z-direction. The dot 
product of the unit vector in the z direction, k, and the 
buoyancy (RaPr) term in Eq. (2) indicates that buoyancy only 
acts in the vertical direction. Therefore the Rayleigh number 
occurs only in the z-momentum equation. 

The relationship showing how viscosity depends on 
temperature, p(T), is given by Eq. (4), This equation is a first 

i 

order expansion of ji(T) expanded about T^ where e = — | 

3T 

Substituting Eq. (4) into the momentum equation, and 
expanding the viscous terms in the momentum equation leads 
to equation (5). 


M-(T) = l + e(T-T b0 ) (4) 

^ = - VP + k • RaPr(T -T b0 ) - Prej^VT • VU + VT • ( VC) T j + 

Pr(l + e(T-T b0 ))v 2 U (5) 


When £ = 0, Eq. (5) reduces to the constant viscosity 
momentum equation. The effect of the parameter, e, on the 
stability of the quiescent base state is the primary focus of this 
paper. 

The nondimensional boundary conditions are given by 
Eqs. ( 6 ) and (7). Equations ( 6 (a) to (c)) represent the no-slip 
conditions and impenetrable wall condition at z = 0. Equation 
( 6 (d)) is the constant temperature condition along the wall. 
The normal stress boundary condition reduces to (Eq. 7(a)) 
when the free surface at z = 1 is assumed to be flat Boundary 
condition (Eq. 7(b)) is the heat flux balance at the free surface, 
where Q* is the dimensional surface heat flux to the 
environment and k^ is the fluid thermal conductivity. Equation 
7(c) is the vector equation for the tangential force balance 
along the free surface and it can be resolved into two equations 
in the x and y directions. 

at z = 0; 0(0) = (U x ,U y ,U z ) = 0; T(0) = T b0 ( 6 (a) to (d)) 


at z = 1; U 2 (1) = 0; 3 -+ 7 T—O, 
az k 0 b 


(7(a),(b)) 


-H(T) 


au 


z + au*^ 


5x a z 


1 + 


au 3u ' 

*- + L 

dy dz 


MaV„ T 


(7(c)) 


The operator, V n , is the surface gradient defined as 

*3 a 9 a a 

i H j — where, i and j are unit vectors in the x and y 

dx dy 

directions respectively. The Marangoni number, Ma, which 

B d Vo 

occurs in Eq. 7(c), is defined as: Ma= — ; — — , where 

k 0 Ho 
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parameter, y' 0 , is defined as 


9a * 


9T* 


, and is often referred to 


as the temperature variation of surface tension (Nield 1964 and 
Adamson 1967) or the differential coefficient of surface 
tension change with temperature (Scriven and Stemling 1964). 
The surface tension, a*, does not appear in our equations or 
boundary conditions since we have assumed a flat interface. 
Further discussion of the nondimensionalization of the free 
surface boundary conditions is found in Scriven and Stemling 
(1964), and Koschmeider (1993). 

The surface heat flux, Q*, has to be expressed in a form 
that is suitable for linearizing the heat flux boundary condition, 
Eq. (7). This is accomplished by expanding Q* about the 
base state surface temperature, T^. The first order expansion 
is given by Eq. (8). As previously noted, the base state varies 
only in the z-direction. Therefore, Q*(T bs ) can be re-expressed 
using Fourier’s law by Eq. (9). 


q*(t*)=q*(t;) + ^- 


( 8 ) 


T* 




q*(t;)=-r; 


. dT 


dz 


= koP 


(9) 


z =d 


Substituting Eq. (8) into Eq. (7(b)), using k^P in place of 

Q*(T bs ) and defining h* = , * the dimensionless heat 

oT T * 

flux boundary condition becomes: 


|Ul + Bi s (T s -Tj = 0 


( 10 ) 


h d 

The dimensionless group, Bi s , is defined as Bi s = — — 

k o 

and is referred to as either the surface Biot number (Pearson, 
1958 and Nield, 1964) or the surface Nusselt number (Scriven 
and Stemling, 1964). 

We note that the three-dimensional mass, momentum, and 
energy equations are given in Eqs. (1) to (3), yet the boundary 
conditions are only specified in the z-direction. After the 
governing equations and boundary conditions are linearized 
and simplified using some vector operations, it is clear that the 
x and y components do not affect the stability of the base state. 
Equations ((1), (3), (5), and (6)), and Eqs. (7(a) and (10)) make 
up the system which we will linearize in the next section. 


Linearize Governing Equations 

The dependent variables are written in terms of the 
following base flow and perturbation variables: 

U = u, T = T b +6, AP = AP b + Ap 

After substituting for T b and VT b , the disturbance equations 
for temperature dependent viscosity become: 


^- = -Vp + k-RaPr0 
ot 


-Pre^k -^Vu + (Vu) T j 

+Pr(l-ez)V 2 u 




(ID 


( 12 ) 


k is the unit vector in the z-direction shown in Fig. 1. The curl 
operator is applied twice to the momentum equation, Eq. (11), 
which yields Eq. (13). 


a(v 2 u) 


9t 


= k 


RaPrV 2 / e-2Pre 


9(v 2 u) 


dz 


+ Pr(l-ez)V 4 u 
(13) 


The first curl operation yields the vorticity equation and 
eliminates the pressure terms. The second curl operation 
decouples the x and y momentum equations from the z- 
momentum and the energy equations. The z-momentum and 
energy equations remain coupled through the buoyancy term in 
Eq. (13), the convective term in Eq. (12), and the tangential 
free surface boundary condition (discussed below). Further- 
more, the relevant stability parameters, Ma and Ra, do not 
appear in either the x or y momentum equation or their 
associated boundary conditions. Given these considerations, 
Eq. (13) reduces to a scalar equation in u z , Eq. (14). 


3(v 2 uJ 


9t 


= RaPr V 2 / 6-2Pre 


a(v 2 uj 


9z 


+ Pr(l-ez)V 4 u 2 


(14) 


The boundary conditions for the perturbed variables 
associated with Eqs. (12) and (14) are given by Eqs. (15) 
and (16). 

9u. 


at z = 0, u z = 0 ; 


9z 


= 0 ; 


0 = 0 


(15(a) to (c)) 


I? 
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atz = 1, U 2 =0; — +Bi s 0(l) = O (16(a) and (b)) 

oz 

-(l-e)^V 2 u z = MaV 2 0 (16(c)) 

Small Disturbance Analysis 

Since Eqs. (12) and (14) are linear, we assume solutions 
for u 2 and <|> are of the form: 

„ v ifa x x + a v y)+Xt ifa x x+a v y]+Xt 

u 2 =w(z)eV x y J ) and 0 = <) > ( z )e ^ x yJ ' 

a x and oty are the dimensionless wavenumbers in the x and y 
directions, and X is the dimensionless frequency. Substituting 
these into Eqs. (12) and (14) results in the following ordinary 
differential equations. 

X<j>(z)-D 2 <[>(z) + (X 2 <J>(z)- w(z) = 0 (17) 

x(D 2 w-a 2 w(z)) = -RaPra 2 <j>-2Pre 

x(D 3 w-a 2 Dw) + Pr(l-ez)(D 4 w-2a 2 D 2 w+a 4 w(z)) 

(18) 

where D =— and a 2 =a 2 +a 2 
dz y 

The boundary conditions at z = 0 reduce to: 

w(0)=0, Dw(0)=0, <t>(0)=0. (19(a) to (c)) 

At z = 1, the flat interface condition, heat flux condition, and 
tangential stress boundary condition are: 

w(l) = 0, D<t>(l) + Bi s <t>(l) =0, (20(a),(b)) 

(l-e)D 2 w = -a 2 Ma<|>(l) (20(c)) 

Equations (17) to (20) are solved to determine whether the 
velocity and temperature disturbances grow or decay for given 
values of the relevant parameters, Ma, Ra, e, and a. This 
problem is also referred to as a temporally developing flow 
problem since the disturbance growth or decay is in time. For 
temporally developing flows, a x and a y are real and the 
eigenvalue, X, is complex. If the real part of X is positive the 
disturbance grows; if the real part of X is negative the 
disturbance decays in time; and if X is zero, the disturbance 
persists unchanged in time. Analytically determined approxi- 
mate solutions may exist, such as the series solutions of Cloot 
and Lebon (1985) and Nield (1964), but we chose a finite 
difference approach because of its ease of implementation. 

Equations (17) to (20) were discretized using a central 
difference scheme. The discretized governing equations, 
Eqs. (21) and (22), were arranged in the form Az = XBz, which 
is the generalized eigenvalue problem. In this analysis. 


coefficients, a through e are functions of z while f, g, and r are 
constants. 

a ( z i > w i — 2 +Mz i )w i ^ 1 +c(z i )w i +d(z i )w i + 1 

+e(z i )w i + 2 +nt> i = ^( fw i-l + 2 w i +f w i + l) 

( 21 ) 

f ^ + **i +f +i+l + w i“ X W (22) 

B is a nonsingular matrix, so it is possible to reduce the system 
to a regular eigenvalue problem of the form Cz = B -1 A = Xz. 
Assuming a flat interface ensures that B is a tridiagonal matrix 
which can efficiently be inverted using a tridiagonal solver. For 
this investigation, the resolution was 50 points across the fluid 
layer in the z-direction. Since there are two differential 
equations the matrix problem has order 100. The nature of the 
boundary conditions permits this to be reduced to 99. 

RESULTS 

Neutral Stability Curves 

The effects of variable viscosity on the neutral stability 
curves for the limiting cases, Marangoni-Benard and Rayleigh- 
Benard problems are shown in Figs. 2 and 3, respectively. 
Results are presented for an insulated free surface, Bi s = 0, 
since this choice of Bi s yielded the smallest values of Ma c and 
Ra c . With the above heat flux restriction the boundary 
conditions corresponding to the shown results are referred to 
as a rigid and conductive lower surface and a free and insulated 
upper surface. In both figures, the loci of points which give 
zero eigenvalues move upwards for increasingly negative 
values of viscosity slopes, e’s, clearly showing the stabilizing 
effect of the temperature-dependent viscosity. In Table 1 we 
have shown the percentage change of the critical Marangoni 
and Rayleigh numbers, relative to the case of constant 
viscosity. The temperature dependent viscosity has a slightly 
greater stabilizing effect on Ma® b for the Marangoni-Benard 
problem, than on Raj 6 for the Rayleigh-Benard problem. The 
superscript "mb” emphasizes that this is the Marangoni-Benard 
problem, therefore the Rayleigh number is zero. The 
superscript "rb" denotes the Rayleigh-Benard problem 
therefore the Marangoni number is zero. This distinction will 
be important when the coupled buoyancy-surface tension 
problem is discussed. 

It is also of interest to know the nature of the new flow 
which develops when the base flow loses stability. Some 
inkling of this can be obtained from the critical value of the 
wavenumber, a. In Figs. 2 and 3 the critical wavenumber 
appears unchanged by the variable viscosity, however Table 2 
shows that the critical wavenumbers, a c 's, decrease slightly as £ 
becomes more negative. 

The constant viscosity curve, e = 0, was first given by 
Pearson for the case of a rigid and conductive lower surface 
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and a free and insulated upper surface. The simpler problem of 
constant viscosity permits a closed form solution of the 
eigenvalue problem and so formulae for the parameter values 
which give zero eigenvalues can be obtained. At smaller values 
of the wavenumber, a, our numerically computed Ma™ b agrees 
with Pearson's exact solution to three significant figures. But as 
the wavenumber is increased, the accuracy of the numerical 
results decreased. The maximum error in Ma™ b with respect to 
Pearson's results is 0.8% which occurs at a = 6. We also 
compared our results with those of Nield for the 
Rayleigh-Benaid problem. In that case, the critical Rayleigh 
number, Ra£\ from the constant viscosity curve in Fig. 3 
agreed with Nield’s tabulated result within 0.1%. 

The Marangoni-Benard problem with a linear temperature 
dependent viscosity profile and a deformable free surface was 
investigated previously by Cloot and Lebon (1985). They 
obtained an analytical solution for that problem in terms of a 
power series expansion. Cloot and Lebon concluded that a 
viscosity profile which decreases linearly with temperature, is 
destabilizing. Their conclusion contradicts the findings of the 
present analysis which shows the above viscosity profile to be 
stabilizing. We believe that this discrepancy occurred for the 
following reason. While Cloot and Lebon permitted the 
kinematic viscosity to vary in the momentum equation, they 
assumed a constant (dynamic) viscosity in the tangential and 
normal stress boundary conditions. Only the tangential 
boundary condition is important here since we assume that the 
interface is flat. The neutral stability results from Cloot and 
Lebon's, Ma<x, are easily corrected for variable viscosity in the 
tangential boundary condition as follows: 
Ma mb (e) = Ma CL * (l-£). After adjusting their solution in 
this manner we find that the results from their power series 
solution agree with the results presented above. Critical 
Marangoni numbers presented in Table 2 differ from the 
corrected Cloot and Lebon results by less than 0.1%. 

Detailed investigations of the variable viscosity Rayleigh- 
Benard problem were performed by Stengel, Oliver, and 
Booker (1982) and White (1988). Stengel et al., considered 
rigid lower and free upper surfaces in addition to other surface 
conditions; however both surfaces were maintained at constant 
temperatures in all cases. White's configuration consisted of 
rigid and conductive (constant temperature) conditions on both 
surfaces. Since our boundary conditions are slightly different, 
the Rayleigh-Benard results shown in Fig. 3 can only be 
compared qualitatively with the Stengel et al. and White 
studies. Stengel et al. showed that when the viscosity decreases 
exponentially with temperature, corresponding to our £<0, and 
when the ratio of (at T^q), to v m i n (at T*) is small to 
moderate, Vm^/v^ < 3000, the effect is that the fluid is more 
stable. On the other hand they found that at larger viscosity 
variation ratios, the effect of the viscosity variation was to 
cause the fluid to go unstable at smaller values of critical 


Rayleigh number. White (1988) performed an extensiveset of 
experiments to study the onset of convection for the Rayleigh- 
Benard problem with temperature-dependent viscosity. He also 
observed that Ra^ increases with an increase in viscosity 
variation for v max /v min < 0(1000) . The results in Fig. 3 agree 
qualitatively with those of Stengel et al. and White. This 
agreement is expected since an exponential viscosity variation 
is well approximated by a linear variation for small £. White 
also examined the effect of variable viscosity on the critical 
wavenumber and his numerical results show that the critical 
wavenumber decreases slightly with increasing viscosity 
variation for v max /v min <55. This is consistent with the 
behavior of a c observed in Table 2. 



a 


Figure 2 Marangoni-Benard Neutral Stability Curves For 
Different Viscosity-Temperature Slope Values, Ra=0 



Figure 3. Rayleigh-Benard Neutral Stability Curves For 
Different Viscosity-Temperature Slope Values, Ma=0 
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Table 1 -Increase of Critical Marangoni Number and 
Critical Rayleigh Number For Marangoni-Benard and 
Rayleigh-Benard Problems 


e 

Ma c (e) - Ma c (0) 

Ra c (e)-Ra f (0) 

Ma c (0) 

Ra c (0) 

0 . 

0 . 

0 . 

-.1 

6.5 

4.9 

-.2 

12.9 

9.8 

-.3 

19.2 

14.6 

-.4 

25.4 

19.3 

-.5 

31.6 

24.0 


Table 2.-Critical Wavenumbers, Marangoni Numbers 
and Rayleigh Numbers For Different Viscosity Slopes 


e 

a n mb (e) 

Ma rab „(e) 

«Ae) 

Ra^e) 

0.0 

1.991 

79.6 

2.085 

668.3 

-0.1 

1.985 

84.8 

2.080 

701.3 

-0.2 

1.974 

89.9 

2.075 

733.8 

-0.3 

1.965 

95.0 

2.070 

765.8 

-0.4 

1.958 

99.9 

2.066 

797.3 

-0.5 

1.951 

104.8 

2.062 

828.4 


Scaling 

In Fig. 4 we show that when the wavenumber is divided by 
its critical value, a c (e), and the Marangoni number is scaled by 
Ma™ b (e), the neutral stability curves for the Marangoni- 
Benard problem, collapse onto a single curve. This figure also 
includes the scaled results for the Rayleigh-Benard problem, 
where the wavenumber and Rayleigh number are divided by 
a c (e) and Ra^ (e) . This means that with a limited amount of 
information, a c (e) and Ma^ b (e), the neutral stability curve for 
a given viscosity slope, £, can be approximated from 
Pearson's exact solution for the constant viscosity neutral 
stability curve. The same is true for the Rayleigh-Benard 
problem, where the neutral stability curve for e^O can be 
obtained approximately by scaling the constant viscosity 
neutral stability curve with known values of a c (e) and Ra * (e) . 
The maximum error that is incurred in this scaling process for 
the Marangoni and Rayleigh curves occurred for 
e = -0.5 and is respectively 3 % and 2 %. 

The dependence of the critical values of Marangoni and 
Rayleigh numbers and wavenumber on the slope of the 
viscosity temperature relation, £, is shown in Fig. 5. This 
figure reveals that the critical Marangoni and critical Rayleigh 
numbers vary linearly with the viscosity slope, E. For the range 
of £ considered, the a c variation with £ can also be 
approximated as linear. The critical numbers, Ma® b (e) and 
Ra^(£), are scaled by the constant viscosity critical numbers, 
Ma™ b (0) and Ra£ b (0). Linear regression results of the 
Figure 5 data are given by Eqs. (23) to (26). Values for a c (£). 


Ma“ b (e), Ra^(e) to three decimal places are given in 
Table 2. 



a c /a t( £ ) 

Figure 4. Scaled Rayleigh-Benard & Marangoni-Benard 
Neutral Stability Curves. 




Figure 5. Critical Marangoni Number, Critical Rayleigh 
Number, and Critical Wavenumbers as a Function of Viscosity 
Temperature Slope 
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Ma™ b (£) , ^ 

c -~ = 1.000-.633e 

Ma” b (0) 

(23) 

Ra ^ = 1.000-. 479 e 

RafCO) 

(24) 

a“ b = 1.990+. 080 e 

(25) 

a* =2.085+.046s 

(26) 


Coupled Buoyancy and Surface Tension Problem 

The stability of a thin layer of fluid with a temperature 
gradient normal to the free surface with combined buoyancy 
and surface tension driving forces and constant viscosity was 
first solved by Nield. In that work it was shown that the 
buoyancy and surface tension forces "reinforced one another” 
with regard to the onset of convection. In other words, Ma c and 
Ra c are approximately inversely proportional for the constant 
viscosity coupled problem. Nield normalized the critical 
Marangoni Number, Ma c , with the critical Marangoni number 
from the Marangoni-Benard problem, Maj? b (0). Similarly, the 
critical Rayleigh number, Ra c was normalized with respect to 
the critical Rayleigh number from the Rayleigh-Benard 
problem, Ra^(O). The interaction between the scaled critical 
values, Ma c and Ra c , for the coupled problem is shown in 
Fig. 6. In our case, we modified Nield’s approach by scaling 
Ma c (e) and Ra c (e) by Ma® b (e) and Raf (£), respectively. The 
variable viscosity results shown in Fig. 6 are in excellent 
agreement with Nield's constant viscosity results. The greatest 
deviation from Nield’s curve is -0.3% and occurs for £ = -0.5 
results. Equation (27) is a curve fit of the constant viscosity 
numerical data (Nield's curve) and is also shown in Fig. 6. 



Figure 6 Scaled Critical Marangoni Numbers and Rayleigh 
Numbers For Different Viscosity Slope Variations 

Ma c (e) 

Ma c mb (e) 


Equations (23) and (24) can be used to compute values of 
Ma® b (£) and Ra*(£). Then using either Fig. 6 or Eq. (27), 
the critical values, Ma c (£) and Ra c (£) for the corresponding 
variable viscosity coupled problem can be determined. 

The results from the previous section are used to 
investigate the stability of a thin layer for two different fluids, a 
large Prandtl number fluid, silicone oil, and a small Prandtl 
number fluid, mercury. Nield’s curve, which is approximated 
by Eq. (27), is used to determine the "critical" temperature 
difference across the fluid layer, AT C , and the dimensionless 
viscosity slope, £, for varying fluid layer thicknesses, d. For a 
given fluid, the appropriate fluid layer thickness required for 
study of the surface tension dominated, buoyancy dominated, 
or coupled instability is also computed from Nield’s curve. The 
critical temperature difference, AT C , is computed by 
substituting Eqs. (23), (24), and the value of £ into Eq. (27) 
and iteratively solving the resulting relationship between, AT C 
and d. After solving for AT C , the curves in Figs. 7(a) and 8(a) 
were obtained by directly substituting AT C , and d into the 
definitions of Ma c /Ma c mb and Ra c /Ra c rb . The critical 
temperature difference and the corresponding slope of 
viscosity versus temperature, £, are plotted as functions of 
fluid depth, d, in Figs. 7(b) and 8(b). 

Figures 7(a) and 8(a) provide a measure of the controlling 
instability mechanism, buoyancy or surface tension, at different 
fluid layer depths, d. The transition from a surface tension 
dominated instability to a buoyancy dominated instability, in 
1 g, occurs at fluid depths of 7.5 mm for the silicone oil and 
8.5 mm for mercury. We (arbitrarily) impose the limits 
Ra c /Ra c rb < .01 to define the Marangoni-Benard instability 
and Ma c /Ma c mb < .01 to define Rayleigh-Benard instability. In 
1 g, the instability can be characterized as the Marangoni- 
Benard for fluid depths less than 1 mm for both silicone oil and 
mercury. For silicone oil and mercury depths greater than 
7.25 cm and 8.15 cm, respectively, the Rayleigh-Benard 
instability occurs in a lg environment For 10"^g, the 
Marangoni instability occurs for d less than 7.5 mm for silicone 
oil and 8 mm for mercury. At 10 g, the Marangoni-Benard 
instability persists beyond the 3.5 cm fluid depth shown in the 
figures. 

Figures 7(b) and 8(b) show that the magnitude of both AT C 
and the dimensionless viscosity slope, £, decrease with 
increasing fluid depth, d. The AT C curves for the reduced 
gravity case, 10"^ g, lie above the lg curves in the above 
figures. The above observations are expected from the trends of 
the critical Marangoni and critical Rayleigh numbers as well 
as from an examination of Eq. (27). Approximating Nield's 
curve by truncating the 1.055 exponent in Eq. (27) to 1, and 
then solving for AT C , we obtain the following approximate 

expression: AT r oc - From this expression it is 

c,d + c 2 gd 3 

readily observed that AT C decreases with an increasing fluid 
depth and gravity. The behavior of -£ is identical to the AT C 
behavior, with magnitude of £ scaled by (Xq and dp*/dT*. 


Ra c (£) 

Ra?fe) I 


(27) 
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Table 3 Thermal Physical Properties of Selected Fluids 
Fluid Properties Taken At 300 K 
Part (a) 


Fluid 

B 

k* 

W/m-K 

Cp* 

J/kg-K 

IB 

IPlgSM 

968 


1507.2 

.096800 


BBSBI 

8.180 

139.3 

.001688 


Part (b) 


Fluid 

e* (dp*/dT*) 
N-s/m 2 -K 

y* (dc*/<n**) 
N/m-K 

V 

1/K 

Silicone OiP 

-6.54e4 

-.00020 

.00018 

Mercury*^ 

-5.35e-6 

-.00006 

.00096 



Layer Depth, d (on) 


Figure 7. Critical Values for Different Fluid Depths 
for Silicone Oil 




Layer Depth, d (cm) 


Figure 8 Critical Values for Different Fluid 
Depths for Mercury 


The viscosity slopes, £, shown in Figs. 7(b) and 8(b) suggest 
that the linear variable viscosity approximation is valid for 
d > 1 mm for silicone oil, and d> 2mm for mercury. The 
variable viscosity effect is negligible, £ < -.01 , for fluid depths 
greater than 9 mm for both the silicone oil and mercury. 

CONCLUDING REMARKS 

Viscosity which is permitted to vary linearly with 
temperature stabilizes the coupled Marangoni-Benard/ 
Rayleigh-Benard problem and both of the limiting cases, 
Rayleigh number equal to zero and Marangoni number equal to 
zero. A steeper (more negative) viscosity slope results in a 
greater stabilizing effect. Cloot and Lebon's conclusion that 
viscosity is destabilizing for the Marangoni-Benard problem 
needs to be adjusted to reflect the variable viscosity in the 
boundary conditions, by simply multiplying their critical 
Marangoni number results by (1-e). The critical parameter 
values for the two limiting cases, Ma® b and Ra^, were found 
to vary linearly with £. Formulae showing the dependence of 
Ma“ b and Ra* behavior on the viscosity slope, £, are given 
by Eqs. (23) and (24). These relations can be used to generate 
the variable viscosity neutral stability curves from the constant 
viscosity neutral stability curves for both the Marangoni- 
Benard and Rayleigh-Benard problems. 

The variable viscosity (linear temperature dependence) 
results for the coupled Marangoni-Benard and Rayleigh- 
Benard problem collapse onto Nield’s curve by appropriately 
scaling Ra c (£) and Ma c (£). The relative importance of 
buoyancy and surface tension to the Benard instability can be 
explicitly shown as a function of fluid layer depths. Specifying 
the fluid depth also fixes AT C through Eq. (27). The parameter 
£, is proportional to the temperature difference across the fluid 
layer, which equals AT C at the bifurcation point. Since AT C is a 
function of fluid depth, the validity of the small £ 
approximation at the bifurcation point that is, the critical 
parameter value where linear stability is lost, can be examined 
for different values of the fluid layer depth. It was determined 
that the linear variable viscosity approximation is valid for 
d > 1 mm for silicone oil, and d > 2 mm for mercury. When the 
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fluid depth exceeds 9 mm for both fluids, the bifurcation point 
can adequately be determined with a constant viscosity 
assumption. 
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